Azimuth correction for data reconstruction in three-dimensional surface-related multiple prediction

ABSTRACT

A best fitting trace in seismic data is determined for a desired trace to be reconstructed. An azimuth correction is calculated for the azimuth difference between the best fitting trace and the desired trace. The azimuth correction is applied to the best fitting trace to reconstruct the desired trace for 3D surface-related multiple prediction.

CROSS-REFERENCES TO RELATED APPLICATIONS

Not Applicable

FEDERALLY SPONSORED RESEARCH OR DEVELOPMENT

Not Applicable

SEQUENCE LISTING, TABLE, OR COMPUTER LISTING

Not Applicable

BACKGROUND OF THE INVENTION

1. Field of the Invention

This invention relates generally to the field of geophysical prospecting. More particularly, the invention relates to the field of surface related multiple prediction in marine seismic surveys.

2. Description of the Related Art

In the oil and gas industry, geophysical prospecting is commonly used to aid in the search for and evaluation of subterranean formations. Geophysical prospecting techniques yield knowledge of the subsurface structure of the earth, which is useful for finding and extracting valuable mineral resources, particularly hydrocarbon deposits such as oil and natural gas. A well-known technique of geophysical prospecting is a seismic survey.

The resulting seismic data obtained in performing a seismic survey is processed to yield information relating to the geologic structure and properties of the subterranean formations in the area being surveyed. The processed seismic data is processed for display and analysis of potential hydrocarbon content of these subterranean formations. The goal of seismic data processing is to extract from the seismic data as much information as possible regarding the subterranean formations in order to adequately image the geologic subsurface. In order to identify locations in the Earth's subsurface where there is a probability for finding petroleum accumulations, large sums of money are expended in gathering, processing, and interpreting seismic data. The process of constructing the reflector surfaces defining the subterranean earth layers of interest from the recorded seismic data provides an image of the earth in depth or time. The image of the structure of the Earth's subsurface is produced in order to enable an interpreter to select locations with the greatest probability of having petroleum accumulations.

In a marine seismic survey, seismic energy sources are used to generate a seismic signal which, after propagating into the earth, is at least partially reflected by subsurface seismic reflectors. Such seismic reflectors typically are interfaces between subterranean formations having different elastic properties, specifically sound wave velocity and rock density, which lead to differences in acoustic impedance at the interfaces. The reflected seismic energy is detected by seismic sensors (also called seismic receivers) and recorded.

The appropriate seismic sources for generating the seismic signal in marine seismic surveys typically include a submerged seismic source towed by a ship and periodically activated to generate an acoustic wavefield. The seismic source generating the wavefield is typically an air gun or a spatially-distributed array of air guns;

The appropriate types of seismic sensors typically include particle velocity sensors (known in the art as geophones) and water pressure sensors (known in the art as hydrophones) mounted within a towed seismic streamer (also know as a seismic cable). Seismic sensors may be deployed by themselves, but are more commonly deployed in sensor arrays within the streamer.

After the reflected wave reaches the seismic sensors, the wave continues to propagate to the water/air interface at the water surface, from which the wave is reflected downwardly, and is again detected by the sensors. The reflected wave continues to propagate and can be reflected upwardly again, by the water bottom or other subterranean formation interfaces. Reflected waves that reflect more than once are termed “multiples” and are typically treated as noise. A particular category of noise comprises multiples that reflect at least once from the water surface and are called surface-related multiples.

Three-dimensional surface-related multiple prediction (3D SRMP), a process to estimate surface-related multiples in seismic data, is a part of 3D surface-related multiple elimination (3D SRME), which strives to attenuate the surface-related multiples by a prediction-subtraction process. The surface-related multiples are first predicted from the seismic data and then the predicted multiples are subtracted from the seismic data to leave a noise-attenuated signal. A first step in this process comprises constructing a multiple contribution gather for a source-receiver trace, which involves the computation of the convolution of pairs of traces over a spatial area. A second step comprises constructing a predicted multiple trace which contains primarily multiple reflections, which involves stacking all the multiple contribution traces belonging to the source-receiver trace. A third step comprises subtracting the predicted multiple traces from the original seismic data.

Thus, the process of predicting 3D surface-related multiples from seismic data, for the combination of a specific source and a specific receiver location, consists of accumulating the results of convolving traces in pairs over a surface area of possible linkage locations. Typically, however, not all pairs of traces to be convolved at all linkage locations are always readily available from seismic data acquisition. This is due to the under-sampling of source and receiver positions, especially in the cross-line direction, that is inherent in marine seismic data acquisition using towed streamers. Furthermore, practical limitations in positioning and navigation systems, as well as wave currents that cause feathering of the streamers and source configuration, prohibit that sources and receivers are positioned exactly where they are desired.

Several methods exist that aim to regularize data in order to reconstruct missing data. These methods may be applied to generate data for 3-D surface-related multiple prediction. These methods known to the industry are based on NMO (Normal Moveout), full DMO (Dip Moveout)—Inverse DMO, AMO (Azimuth Moveout), or migration operators. These methods have in common that they aim to generate new data from the acquired data that are more similar to the desired traces.

Thus, a need exists for a method for data reconstruction to generate desired traces from the set of acquired traces, in order to achieve more accurate 3D surface-related multiple prediction.

BRIEF SUMMARY OF THE INVENTION

The invention is a method of azimuth correction for data reconstruction in 3D surface-related multiple prediction. A best fitting trace in seismic data is determined for a desired trace to be reconstructed. An azimuth correction is calculated for the azimuth difference between the best fitting trace and the desired trace. The azimuth correction is applied to the best fitting trace to reconstruct the desired trace for 3D surface-related multiple prediction.

BRIEF DESCRIPTION OF THE DRAWINGS

The invention and its advantages may be more easily understood by reference to the following detailed description and the attached drawings, in which:

FIG. 1 is a schematic diagram of traces for data reconstruction for 3D surface-related multiple prediction, as in the method of the invention;

FIG. 2 is a flowchart illustrating an embodiment of the invention for azimuth correction for data reconstruction for 3D surface-related multiple prediction; and

FIG. 3 is a flowchart illustrating another embodiment of the invention for azimuth correction for data reconstruction for 3D surface-related multiple prediction.

While the invention will be described in connection with its preferred embodiments, it will be understood that the invention is not limited to these. On the contrary, the invention is intended to cover all alternatives, modifications, and equivalents that may be included within the scope of the invention, as defined by the appended claims.

DETAILED DESCRIPTION OF THE INVENTION

In general, the invention is a method for transforming seismic data into noise-attenuated data. In particular, the invention is a method for seismic data reconstruction using azimuth correction for surface-related multiple prediction schemes. The invention is used in the prediction phase of a 3D prediction-subtraction process to attenuate surface-related multiples in seismic data. The process of predicting 3D surface-related multiples from seismic data comprises accumulating the results of convolving pairs of traces over a surface area of possible linkage locations for each combination of a selected source location and a selected receiver location.

FIG. 1 is a schematic diagram of traces for data reconstruction for 3D surface-related multiple prediction, as in the method of the invention. The location at which surface-related multiples are to be predicted (and subtracted) is determined by a source-receiver trace 10. A source-receiver trace 10 is given by a combination of a selected source location 11 and a selected receiver location 12. The process of predicting the surface related multiples from seismic data for this source-receiver trace 10 comprises convolving pairs of traces and then accumulating the convolution results over a surface area 13. The first trace in each pair of convolved traces is called a source-side linkage trace 14 and has a source location at the selected source location 11 of the source-receiver trace 10 and has a receiver location at a selected linkage location 15. The second trace 16 in the pair is called a receiver-side linkage trace 16 and has its source location at the same selected linkage location 15 and has a receiver location at the selected receiver location 12 of the source-receiver trace 10.

The linkage locations 15 for all pairs of source-side linkage traces 14 and receiver-side linkage traces 16 are distributed over the surface area 13. Each linkage location 15 corresponds to a potential surface reflection point of the surface related multiples associated with this source-receiver trace 10. The traces created by the convolution of the pairs of source-side linkage traces 14 and receiver-side linkage traces 16 are referred to as multiple contribution traces and the ensemble of these multiple contribution traces (before accumulation) is referred to as a multiple contribution gather. The accumulation of the multiple contribution traces into the predicted surface related multiples comprises a process commonly known as Fresnel stacking. The extent or coverage in the surface area 13 of the surface reflection points as incorporated (after accumulation) in the prediction is commonly referred to as the aperture of the prediction.

The aperture can be chosen as the simplest case, a rectangular surface area 13 surrounding the midpoint 17 of the source-receiver trace 10 for the specific source location 11 and specific receiver location 12. The rectangular surface area 13 can be parameterized by an inline direction 18 and a cross-line direction 19. However, the surface area 13 of the aperture does not have to be a rectangular area or oriented in the inline 18 and cross-line directions 19. This is not a limitation of the invention, but merely employed for simplicity of illustration.

However, as discussed above in the Background section, not all pairs of traces to be convolved at all linkage locations are always readily available from marine seismic data acquisition. The invention is a method for seismic data reconstruction by means of azimuth correction to solve this problem.

A conventional approach is to reconstruct the desired traces from the nearest traces in the acquired seismic data. The seismic data is expressed as (x^(CMP),y^(CMP),h,φ,t), where the acquired seismic data is transformed to this format, if necessary. CMP (Common Midpoint) defines the x and y positions of the location mid-way between the source and receiver locations for the trace. Here, x and y are orthogonal coordinates in a processing grid, typically oriented so that x is in the inline direction and y is in the cross-line direction of the seismic survey, although this orientation can be user-defined. Azimuth, φ, is here defined as the angle between the direction of data acquisition (typically, sail line direction in marine data acquisition) and the straight line between the source and receiver locations. Offset, h, is defined as the total distance between the source and receiver location for the individual trace. Arrival time, t, is the time an event is recorded on the trace.

In this method, the desired trace is selected as the best available trace that minimizes the difference between the best fitting (nearest) and the desired trace. The difference Φ is a weighted sum of the differences in CMP x and y positions, offset h, and azimuth φ between the nearest and desired traces, which can be expressed as Φ=α|x _(b) ^(CMP) −x _(d) ^(CMP) |+α|y _(b) ^(CMP) −y _(d) ^(CMP) |+β|h _(b) −h _(d)|+ε|φ_(b)−φ_(d)|,  (1) where α, β and ε are user-defined weights, and subscripts b and d denote best fitting and desired traces, respectively.

Unfortunately, for the reasons stated above, the desired trace with the correct CMP x and y positions, offset h, and azimuth φ is often not present in the set of acquired traces. Conventional data acquisition often aims at a distribution of sources and receivers such that at each CMP x and y, and offset h, one trace is available. In this case, only one azimuth φ is available for the given CMP x and y, and offset h combination. Further, certain modern data acquisition geometries aim at having more than one trace available at such CMP x, y, and offset h location. In these cases, more than one trace is available, each with their own azimuth φ, for each CMP x and y, and offset h combination.

This means that in conventional practice, the best fitting trace is found by minimizing the difference Φ as a weighted sum of differences in CMP x and y, and offset h only, thereby ignoring azimuth φ. Thus, Equation (1) is applied as the shortened version: Φ=α|x _(b) ^(CMP) −x _(d) ^(CMP) |+α|y _(b) ^(CMP) −y _(d) ^(CMP) |+β|h _(b) −h _(d)|,  (2) without the final difference term in azimuth φ.

After finding the best fitting, or nearest, trace, a differential move-out is applied to this data trace in order to correct for the offset difference between the desired trace and the best fitting trace. This differential move-out operator only corrects for offset differences; it cannot correct for the difference in azimuth between the desired trace and the best fitting trace.

The invention is a method to incorporate a correction for azimuth differences between the desired trace and the best fitting trace in the data reconstruction scheme described above. The best fitting trace is found, an azimuth correction is computed for the desired and best fitting traces, and the azimuth correction is then applied to the best fitting trace at each sample. Because the correction is determined and then applied per trace, and per sample, the method of the invention becomes computationally very efficient. The method of the invention is thus fundamentally different from the data regularization approach described in the Background section. The method of the invention is utilized to apply an azimuth correction during data reconstruction for 3-D surface-related multiple prediction:

FIGS. 2-3 show flowcharts illustrating embodiments of the invention for azimuth correction for data reconstruction for 3D surface-related multiple prediction. FIG. 2 is a flowchart illustrating a general embodiment of the invention. FIG. 3 is a flowchart illustrating a more detailed embodiment of the invention

FIG. 2 is a flowchart illustrating an embodiment of the invention for azimuth correction for data reconstruction for 3D surface-related multiple prediction.

At block 20, a best fitting trace is determined for a desired trace to be reconstructed.

At block 21, an azimuth correction is calculated for the azimuth difference between the best fitting trace and the desired trace.

At block 22, the azimuth correction is applied to the best fitting trace to reconstruct the desired trace for 3D surface-related multiple prediction.

FIG. 3 is a flowchart illustrating another embodiment of the invention for azimuth correction for data reconstruction for 3D surface-related multiple prediction. FIG. 3 is a more general embodiment of the method of the invention discussed with reference to the flowchart in FIG. 2.

At block 30, marine seismic data is obtained in a (x^(CMP),y^(CMP),h,φ,t) format, where the acquired seismic data is transformed to this format, if necessary.

At block 31, for each data sample from block 30, with CMP coordinates x^(CMP) and y^(CMP) and time sample t, the predominant dip in the x and y directions is determined. The dips can be determined by any appropriate method known in the art. The dips can be determined from any appropriate subset of the total acquired data volume, such as, for example, from common offset volumes. These determinations generate a dip-estimate volume.

Dip is the slope between two adjacent traces in the gather of the same seismic event. The dip is measured as a ratio between a time difference and a space difference (analogous to dt/dx). The dip is typically decomposed into two components, oriented along the x and y axes of the processing grid being used. Although a dipping event is generally meant to be non-horizontal, a horizontal event is not excluded.

The method of the invention is not restricted by a need to estimate dips for each desired trace to be reconstructed. In an alternative embodiment, dip estimation may be applied on a sparse grid where dip estimates for intermediate positions are obtained by interpolating dip estimates acquired on the sparse grid.

In another alternative embodiment, a finite number of predominant dips may be estimated for any data sample. This multiple dip estimation can be achieved by separating different dip components prior to dip estimation using well-known methods like f-k (frequency-wavenumber) or Radon filtering, or by utilizing well-known dip-decomposition methods, such as slant-stack filtering. The final azimuth correction to be applied will then be a weighted sum of the individual corrections for each dip.

At block 32, a filter, for example a median filter, is applied to the dip-estimate volume from block 31 to smooth discontinuities, if desired. This step is optional.

At block 33, a set of traces to reconstruct is selected. This set of traces comprises the traces that require reconstruction to obtain desired traces for 3-D surface-related multiple prediction.

At block 34, for each of the traces in the set of traces from block 33, a best fitting trace from the seismic data from block 30 is determined. In one embodiment, this best fitting trace is determined by application of Equation (2). However, this is not a restriction of the invention. In another embodiment, a trace from the seismic data from block 30 is determined that is close to the trace to be reconstructed, but not necessarily a best fitting trace until after the azimuth correction is applied in block 36 below.

At block 35, an azimuth correction is determined for the desired trace from block 33 and the best fitting trace from block 34. The azimuth correction represents the difference between the azimuth of the best fitting trace from block 34 and that of the desired trace to be reconstructed from block 33. The azimuth correction is applied to each time sample t in the best fitting trace from block 34. As discussed above, the azimuth correction can be a weighted sum of azimuth corrections for each dip at the data sample.

In one embodiment, the azimuth correction is determined based upon half-offset, arrival time, slowness, dip, and dip direction. In one embodiment, the azimuth correction φ_(cor) used is given by:

$\begin{matrix} {{\phi_{cor} = {\frac{\left( h^{\prime} \right)^{2}}{2\; t}{p^{2}}{\sin\left( {\phi_{b} - \phi_{d}} \right)}{\sin\left( {\phi_{b} + \phi_{d} - {2\; p_{\theta}}} \right)}}},} & (3) \end{matrix}$ with magnitude of slowness p given by:

$\begin{matrix} {{p} = \sqrt{\left( \frac{\mathbb{d}t}{\mathbb{d}x} \right)^{2} + \left( \frac{\mathbb{d}t}{\mathbb{d}y} \right)^{2}}} & (4) \end{matrix}$ and slowness direction p_(θ) given by:

$\begin{matrix} {{p_{\theta} = {\cos^{- 1}\left( \frac{\frac{\mathbb{d}t}{\mathbb{d}x}}{p} \right)}},} & (5) \end{matrix}$ where h′ is the half-offset, t is the arrival time, p is slowness, and θ is the dip direction.

At block 36, the azimuth correction from block 35 is applied to the best fitting trace from block 34 to reconstruct a corresponding desired trace for 3D surface-related multiple prediction. Because the azimuth correction factor is applied per trace and per sample, after the selection of the best fitting trace, the method is computationally very efficient.

In addition, the azimuth correction as determined by the method of the invention is independent of wavefield propagation velocities, Thus, no a priori knowledge of wavefield propagation velocities is required in Equation (3) above. However, this velocity-independent approach is a possible advantage, but not a required limitation of the invention. In other embodiments, other equations could be derived that achieve the same result, an azimuth correction, but would require more detailed subsurface information, such as wavefield propagation velocities. For example, the following equation could be applied to determine a time correction Δt_(h):

$\begin{matrix} {{{\Delta\; t_{h}} = {\frac{2}{v^{2}t_{1}}2\; h_{1}\Delta\; h\;{\sin^{2}(\psi)}{\cos^{2}\left( {\phi_{1} - \theta} \right)}}},} & (6) \end{matrix}$ where v is velocity of the medium, t₁ is a first traveltime, Δh=h₂−h₁ is offset difference for first offset h₁ and second offset h₂, ψ is dip angle, φ₁ is a first source-receiver azimuth, and θ is dip direction. In another embodiment, the following equation could be applied to determine a time correction Δt_(φ):

$\begin{matrix} {{{\Delta\; t_{\phi}} = {\frac{2h^{2}}{v^{2}t_{1}}{\sin^{2}(\psi)}{\sin\left( \phi_{2} \right)}{\cos\left( {\phi_{2} - {2\theta}} \right)}}},} & (7) \end{matrix}$ where φ₂ is a second source-receiver azimuth.

In another embodiment, the method of the invention can also be utilized in the process of predicting 3-D interbed multiples rather than surface-related multiples. In 3D IMP, (Interbed Multiple Prediction), a triplet of traces and two convolutions are used as part of the method, instead of two traces and a single convolution as in 3D SRMP. The method of the invention can be used to reconstruct the missing traces for each of the two convolutions.

The invention has been discussed above as a method, for illustrative purposes only, but can also be implemented as a system. The system of the invention is preferably implemented by means of computers, in particular digital computers, along with other conventional data processing equipment. Such data processing equipment, well known in the art, will comprise any appropriate combination or network of computer processing equipment, including, but not be limited to, hardware (processors, temporary and permanent storage devices, and any other appropriate computer processing equipment), software (operating systems, application programs, mathematics program libraries, and any other appropriate software), connections (electrical, optical, wireless, or otherwise), and peripherals (input and output devices such as keyboards, pointing devices, and scanners; display devices such as monitors and printers; computer readable storage media such as tapes, disks, and hard drives, and any other appropriate equipment).

In another embodiment, the invention could be implemented as the method described above, specifically carried out using a programmable computer to perform the method. In another embodiment, the invention could be implemented as a computer program stored in a computer readable medium, with the program having logic operable to cause a programmable computer to perform the method described above. In another embodiment, the invention could be implemented as a computer readable medium with a computer program stored on the medium, such that the program has logic operable to cause a programmable computer to perform the method described above.

It should be understood that the preceding is merely a detailed description of specific embodiments of this invention and that numerous changes, modifications, and alternatives to the disclosed embodiments can be made in accordance with the disclosure here without departing from the scope of the invention. The preceding description, therefore, is not meant to limit the scope of the invention. Rather, the scope of the invention is to be determined only by the appended claims and their equivalents. 

We claim:
 1. A method for geophysical prospecting, comprising: disposing a marine seismic streamer in a body of water; responsive to signals recorded at receivers in the marine seismic streamer, calculating seismic traces representing physical pressure wavefields traveling in the body of water incident on the receivers; correcting for irregular source and receiver positioning by transforming the seismic traces to produce reconstructed desired traces representing pressure wavefields as they would appear if incident on regularly positioned receivers from regularly positioned sources, the transforming comprising: determining a best fitting trace in the seismic data for a desired trace to be reconstructed; calculating an azimuth correction for the azimuth difference between the best fitting trace and the desired trace; applying the azimuth correction to the best fitting trace to reconstruct the desired trace for 3D surface-related multiple prediction; and recording the desired traces.
 2. The method of claim 1, wherein the calculating an azimuth correction comprises applying the following equation: ${{\Delta\; t_{h}} = {\frac{2}{v^{2}t_{1}}2\; h_{1}\Delta\; h\;{\sin^{2}(\phi)}{\cos^{2}\left( {\psi_{1} - \theta} \right)}}},$ where v is velocity of the medium, t₁ is a first traveltime, Δh=h₂−h₁ is offset difference for first offset h₁ and second offset h₂, φ is dip angle, ψ₁ is a first source-receiver azimuth, and θ is dip direction.
 3. The method of claim 1, wherein the calculating an azimuth correction comprises applying the following equation: ${{\Delta\; t_{\psi}} = {\frac{2h^{2}}{v^{2}t_{1}}{\sin^{2}(\phi)}{\sin\left( \psi_{2} \right)}{\cos\left( {\psi_{2} - {2\theta}} \right)}}},$ where h id offset, v is velocity of the medium, t₁ is a first traveltime, φ is dip angle, ψ₂ is a second source-receiver azimuth, and θ is dip direction.
 4. The method of claim 1, wherein the calculating an azimuth correction comprises: calculating an azimuth correction for the best fitting trace and the corresponding desired trace, based upon half-offset, arrival time, slowness, dip, and dip direction.
 5. The method of claim 4, wherein the calculating an azimuth correction comprises applying the following equation: ${\phi_{cor} = {\frac{\left( h_{2} \right)^{2}}{2\; t}{p^{2}}{\sin\left( {\phi_{b} - \phi_{d}} \right)}{\sin\left( {\phi_{b} + \phi_{d} - {2\; p_{\theta}}} \right)}}},{with}$ ${p} = \sqrt{\left( \frac{\mathbb{d}t}{\mathbb{d}x} \right)^{2} + \left( \frac{\mathbb{d}t}{\mathbb{d}y} \right)^{2}}$ and ${p_{\theta} = {\cos^{- 1}\left( \frac{\frac{\mathbb{d}t}{\mathbb{d}x}}{p} \right)}},$ where h₂ is the half-offset, t is the arrival time, p is slowness, φ_(b) and φ_(d) are best fit azimuth and desired azimuth, respectively, and θ is dip direction.
 6. The method of claim 1, wherein the determining a best fitting trace further comprises: obtaining seismic data in a (x^(CMP),y^(CMP),h,φ,t) format; determining a predominant dip in x and y directions for CMP coordinates x^(CMP) and y^(CMP) and arrival time tin each trace in the seismic data; applying a filter to the dips to smooth discontinuities; selecting a set of desired traces to reconstruct; and determining a best fitting trace from the seismic data for each desired trace in the selected set of traces.
 7. The method of claim 6, wherein the applied filter is a median filter.
 8. The method of claim 6, wherein the determining a best fitting trace comprises minimizing a difference Φ in a weighted sum of differences in CMP x and y positions x^(CMP) and y^(CMP) offset h, and azimuth φ between the best fitting and desired traces.
 9. The method of claim 8, wherein the minimizing a difference Φ comprises applying the following equation: Φ=α|x _(b) ^(CMP) −x _(d) ^(CMP) |+α|y _(b) ^(CMP) −y _(d) ^(CMP) |+β|h _(b) −h _(d)|+ε|φ_(b)−φ_(d)|, where α, β and ε are user-defined weights, and subscripts b and d denote best fitting and desired traces, respectively.
 10. The method of claim 8, wherein the minimizing a difference Φ comprises applying the following equation: Φ=α|x _(b) ^(CMP) −x _(d) ^(CMP) |+α|y _(b) ^(CMP) −y _(d) ^(CMP) |+β|h _(b) −h _(d)|  (1) where α and β are user-defined weights, and subscripts b and d denote best fitting and desired traces, respectively. 